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ABSTRACT 

This report describes a simplified method of computing Lagrangian coefficients 
for curvilinear interpolation, which may be used when tables of Lagrangian 
coefficients are not available or when tables are available but the coefficients are 
not tabulated for the exact fraction of the interval to which the interpolation is 
to be made. 
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I. INTRODUCTION 

The common formulas of curvilinear interpolation, 3 such as those 
of Newton, Gauss, Stirling, Bessel, and Lagrange 4 when carried to 
completion, are all equivalent in that they approximate the function 
by a polynomial of suitable degree (the number of points required is 
one greater than the degree of the polynomial assumed). All of these 
formulas except that of Lagrange require, as a first step, the compu- 
tation of a "table of differences" for the function, to an order equal 
to the degree of the polynomial. Lagrangian interpolation, on the 
other hand, is carried out by multiplying each of the tabulated values 
on which the interpolation is to be based by a suitable coefficient — 
the Lagrangian coefficients — and taking the sum of these products. 
This "cumulative" multiplication can be carried out as a continuous 
operation on a modern calculating machine, and when such a machine 
is available Lagrangian interpolation is usually more rapid than the 
other methods referred to. 

Very extensive tables of Lagrangian interpolation coefficients are 
now available 5 (references to other, less complete, tables are given in 
this work). 

1 This investigation was performed at the National Bureau of Standards as part of the work of the Ameri- 
can Petroleum Institute Research Project 44 on the "Collection and Analysis of Data on the Properties of 
Hydrocarbons." 

2 Research Associate on the American Petroleum Institute Research Project 44 at the National Bureau 
of Standards. 

8 The discussion in this report is restricted to interpolation in tables in which a function is tabulated at 
equally spaced intervals of the independent variable. 

4 David, Gibb, Edinburgh Mathematical Tracts No. 2, Interpolation and numerical integration (G. 
Bell & Sons, Ltd., London, 1915). 

• Tables of Lagrangian interpolation coefficients, Mathematical Tables Project, Works Project Admin* 
istration, Federal Works Agency (National Bureau of Standards; Columbia University Press, New York, 
N. Y., 1944). 
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In this report there is described a simplified method of computing 
the Lagrangian coefficients, which may be used when-tables of La- 
grangian coefficients are not available, or when tables are avail- 
able, but the coefficients are not tabulated for the exact fraction of 
the interval to which the interpolation is to be made. This method 
requires only a single division for the computation of each coefficient, 
for interpolation on the basis of any number of tabulated values. 
The coefficients, as thus computed, are not normalized — they may be 
normalized in the usual way, although it is usually more convenient 
to interpolate by means of the non-normalized Lagrangian coefficients. 

Lagrangian interpolation, with the use of tables of Lagrangian 
coefficients (see footnote 5), or by the method described in this report, 
is recommended as a rapid and generally useful method of interpola- 
tion in those tables in which a function, or property, is tabulated at 
several equally spaced values^fofjthe independent variable. The 
number of tabulated values, or the ^degree of the polynomial, on 
which the interpolation is based must be chosen to represent the 
function with sufficient accuracy in the interval in question. 

II. GENERAL METHOD 

The form of the modified Lagrangian interpolation formula will be 
stated, and then the derivation of the formula will be given. 

Let the (n+1) tabulated values of the function y(x) on which the 
interpolation is to be based be y , y u . . ., y k , . . ., y nj and let the 
corresponding values of the independent variable — equally spaced at 
intervals of h — be x , x ir . . , x k , . . . , x n . The (n+ 1) values of y 
determine a polynomial of the ftth degree in x, P(x), which is an approx- 
imation to the function y(x) in the range x to x n . The value of P(x), 
and therefore, approximately, of ?/, at a given value of x y is given by 
the modified Lagrangian formula 

P(x)=i>*y*/i>*, (i) 

k~Q *«0 

where the coefficients, a k , are given by 

a*=±(j)/|fe-x)/A|. (2) 

\(x k — x) /h\ is the absolute value of (x k —x)lh, and (J) is the binomial 
coefficient 

(j)-n!/*!(n~fc)L (3) 



The signs of the coefficients, a k , are determined as follows: If x 



lies between x m and x m+1 
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then the coefficients, . . . a m _ 2 , a m , a m+u a m+3 , . . ., are positive, and 
the coefficients, . . ., a m _ 3 , a m _i, a OT+2 , &m-H, • • •> are negative. The 
pattern of signs is thus of the type 

...- + - + - + + -+ + 
± 



POINT OF INTERPOLATION 



The coefficients, a ky are proportional to the Lagrangian interpolation 
coefficients but are not normalized. In fact, the sum of the a k s is 
given by the expression 

2>*=*U/ n \(x k -x)/h\, (4) 

which may be used as a check ou the calculation of the a k s. 

If several interpolations to the same fraction of an interval are to 
be made, it is more convenient to put eq 1 in the form 

P(x)=itA k y k , (5) 

where the ^IVs are the normalized Lagrangian coefficients, 

n 

A k =a k /Jla k , (6) 

which have been tabulated elsewhere (see footnote 5). As used 
here, "normalized" indicates that the algebraic sum of the coefficients 
is unity. 

These results may be derived in the following way: P(x) may be 
written in the form 

p ( T \ __. y^l \ x x o)\ x x v - • • \ x x k -i)(x x k +i) . . . (x x n ) | 
£oLfe— x o)( x k— x l) . . . Oa— ^-i)fe— x k+1 ) . . . {x k — x n )J y * 

(7) 

as P(x) (by definition) and the expression on the right of eq 7 are both 
polynomials of the nth degree in x } and are equal to y k when x=x k , for 
k=0, 1, . . . , n. But two polynomials of the nth degree which 
coincide at (n-\-l) points must be identical, which establishes eq 7. 
Rearrangement of eq 7 yields 

KnPXnrD • • • (?)/-] SE-AV'XV) ■ ■ • 
\ h A a )• • ■\Jr)] y "i\nr) w 
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The expression in brackets under the summation sign is equal in 
absolute value to the binomial coefficient, (J). Therefore if the signs 
of the terms are temporarily disregarded, as indicated by the ± sign, 
eq 8 may be rewritten as 

P{x) =[i U |(*-z*)/a| ]§[±(*)/|<**-*)A|]v*- W 

Inspection of eq 7 shows that the signs of the terms follows the rule 
already given for the signs of the coefficients a k . In view of the defini- 
tion of a k , eq 2, therefore, eq 9 becomes 



tl n I "1 n 

—. n (x—x k )lh\ 5>a#*. 



(10) 



In the special case in which |fe=l, for Jc=0, 1, . . . , n, the polynomial 
P(x) must have the value unity for all values of x. It follows that 



r i n n * 

- n (x-x k )lh\ 2>*= 



(ii) 



Equations 10 and 11 lead immediately to eq 1, and eq 11 is equivalent 
to eq 4. 

III. NUMERICAL EXAMPLE 

Given the following (equally spaced) values of the descending 
exponential function, e~ x , find by interpolation the value of e~ x at 
z=0.54316: 



X 


y=e-* 


z = 0.52 
Z! = 0.53 
z 2 = 0.54 
z 3 = 0.55 
z 4 = 0.56 


2/0=0.594 520 548 
2/1 = 0.588 604 969 7 
2/ 2 = 0.582 748 252 4 
2/3 = 0.576 949 810 4 
2/4 = 0.571209 063 8 



Solution: The interval, A, of the independent variable, x, is 0.01. 
The values of the binomial coefficients for n=4 are 1, 4, 6, 4, and 1. 
The value x=0. 54316 lies between x 2 =0.54 and £3=0.55. Then the 
values of the a k s, or the non-normalized Lagrangian coefficients, are 
given by eq 2: 

ao= + (£)l\(xo—x)/h\ = + 1/2.316 = + 0.431778 9 
a 1 =^-({)/\(x 1 -x)/h\ = -4/1.316- - 3.039 513 7 
a 2 = + (t)l\(x 2 —x)/h\ = +6/0.316- + 18.987 341 8 
a 3 = + (S)/| (x 3 — ») A| = +4/0.684- + 5.847 953 2 
a 4 =-(t)/\(x 4 -x)/h\ = - 1/1.684- - 0.593 824 2 



23 ^=21.633 736 



mmw^mmmmmmmm 
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A check may be obtained from eq 4: 



/ n 



(x k —x)/h 



=24/(2.316) (1.316) (0.316) (0.684) (1.684) 
=21.633 736 0. 



The value of P(x) at £=0. 54316 may now be calculated from eq 1 
(the cumulative multiplication involved may be performed as a 
continuous operation on the calculating machine) : 

4 / 4 

£=0 / &=0 

/ 12.567 246 5 X 
V21.633 736 07 

= 0.580 909 673. 

The true value of y=e~ x , at x= 0.54316, correct to 10 decimal 
places, is 

?/=0.580 909 674 4. 

When the same coefficients are to be used several times, it is con- 
venient to compute the normalized Lagrangian coefficients, or the 
A k s, from eq 6: 

A = +0.019 958 592 

A x =— .140 498 789 

A 2 =+ .877 672 807 

Az=+ .270 316 380 

A±=— .027 448 990 



ZM*= 1.000 000 000 



These coefficients differ by one or two units in. the last place from the 
correct values (see p. 232 of the reference in footnote 5) because the 
a k s were computed to only seven decimal places. The value of P(x) 
at #=0.54316 may now be calculated from eq 5: 

P(x)^A k y k 

= 0.580 909 674. 
Washington, March 3, 1945. 

o 



